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Abstract 



We consider a bistable reaction-diffusion system on the interface of deep 
fluid interacting with Marangoni flow. The method of matched asymptotic 
expansions is used to resolve the singularity at a sharp interface between the 
alternative states, and to compute the self- induced flow velocity advecting 
the domain boundary. It is shown that Marangoni flow serves as an efficient 
mechanism preventing the spread of the state with a higher surface tension 
when it is dynamically favorable. 
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It has been long realized that chemical instabilities can both induce and be strongly af- 
fected by hydrodynamic flows. In early experiments with BZ and other oscillatory reactions 
this interaction led to a transition from propagating waves in confined systems to cellular 
structures in systems with a free interface where Marangoni flow was excited by concenra- 
tion gradients of chemical origin Q. Strong interaction between chemical and interfacial 
instabilities has been the subject of a number of experimental studies 0. A related direc- 
tion was the study of chemical fronts affected by buoyancy |§ ; both effects were apparently 
relevant in observations of patterns of photochemical origin |4j]. 

Following the linear analysis of chemical instability coupled with Marangoni flow ||, 
theoretical analysis of the nonlinear development of reaction-convection instabilities has 
been restricted so far either to very shallow layers enabling lubrication approximation ||, 
or to systems close to a bifurcation point ; both approaches are inapplicable to the most 
practical case of sharp reaction fronts. The only other alternative was numerical modeling 

§■ 

The most powerful technique for constructing and analysing non-equilibrium patterns 
of chemical origin is based on tracing the motion of sharp interfaces between alternative 
states of a bistable system interacting with a long-range field |gjn|. Maran Efoni flow due 



to concentration gradients in the transitional region separating the two alternative states 
and, in its turn, advecting this boundary, well qualifies as a long-range controlling agent; its 
action is, however, highly nonlocal, which poses serious analytical difficulties. 

Recently, boundary integral technique was applied to computation of the motion of 



interfacial domains coupled to creeping motion of the underlying fluid [|iT[|12| . The effect of 
a sharp interface was, however, the emergence of a finite time singularity Jl2[ that ruled out 
formation of stable structures. I shall show in this Letter that the singularity can be indeed 
resolved when finite, though low diffusivity is taken into account. The method of matched 
asymptotic expansions (akin to that used for the analysis of motion of vortex lines [H|) will 



be applied to resolve the structure of the domain boundary, and to obtain a finite stationary 



velocity of the induced flow. I shall further show that inward flow can stabilize a solitary 
spot of a dynamically prevalent state with a higher surface tension. 

Consider creeping motion of a viscous incompressible fluid in an infinitely deep and 
infinitely extended layer z < induced by the Marangoni force due to a given stationary 
distribution of surface tension er(x) on the undeformable free boundary z = 0. Since the 
vertical vorticity vanishes, the 2d horizontal velocity vector v and the vertical velocity w 
are expressed, respectively, as v = V^z, w = — V 2 x, where V is the 2d vector differential 
operator. The velocity potential x( x , z) is determined by solving the Stokes equation 

(v 2 + a 2 ) 2 x = o (i) 

with the Marangoni boundary conditions at the free surface z = 0: 

X (x,0) = 0, f i Xzz (x,0) = -V(T, (2) 

where \x is dynamic viscosity. 

The solution is found most readily by Fourier transform yielding 

m = -^ zekz& ( k )> ( 3 ) 

where hats denote Fourier transforms of the respective functions, and k — |k|. The only 
quantity relevant for our purpose is the interfacial velocity, which is computed as v = V0, 
where the interfacial flow potential 0(x) = Xz( x ) 0) can be expressed in the form 

<t>(x)= t T 1 J G(|x-flM£)de, (4) 

where the integration is carried out over the entire free surface. The kernel G(r) depends 
only on the 2d distance r = |x — £|, and is obtained as the inverse transform of x 2 (0): 

G(r) = —= / fe- 1 e- ik - x rf 2 x = — / J (kr)dk = — . (5) 
v ; 8vr 2 J 4tt Jo uv ' Anr v ' 

It is convenient to express the velocity through the gradient of surface tension: 

v(x)= M - 1 | G(x-0W(£K. (6) 



The solution exists provided cr(x) = const at |x| — > oo. 

The distribution of the surface concentration 9 of an insoluble surfactant on the free sur- 
face obeys the convective reaction-diffusion equation, which we shall write in a dimensionless 
form 

e t = \v 2 e- v-(c0) + f{9). (7) 

Here f(9) is the dimensionless net surfactant source due to chemical reactions and exchange 
with the gas phase; the time scale is the characteristic reaction time r; the coordinates are 
scaled by the diffusional length y/2Dr, where D is the surface diffusivity; c = vyV/2.0 is 
the dimensionless velocity. We are interested in the case when at least two stable stationary 
solutions exist; therefore the function f(9) should have three zeroes, and its derivative at the 
smallest and largest zero should be negative. The simplest function possessing this property 
is a cubic. By shifting and rescaling, the two stable zeroes can be placed at 9 = and 9 = 1. 
Thus, we can assume f{9) = —9(9 — q){9 — 1), where < q < 1. 

If c = const, Eq. (0) with the cubic source has a simple stationary solution 

9 = = - ( 1 + tanh -) , (8) 

corresponding to the velocity c = — cq = \ — q directed along the x axis. The front translates 
normally to itself towards the upper state with the "chemical" speed Co in the coordinate 
frame comoving with the local interface flow velocity c. Hence, cq — c ■ n is the displacement 
speed of a point on the front along the normal n. This speed can be corrected to account for 
the curvature effect. The front is stationary when the interface is advected in the direction 
opposite to chemical propagation with the same speed. 

In the presence of Marangoni flow, Eq. (0) becomes an integro-differential equation, 
where c is expressed with the help of Eq. (H) as 

c(x) = -M|G(x-OV0(£)de, (9) 

where M = (Act / r /2D) is the Marangoni number based on the characteristic diffusional 
length; Act is the difference of surface tension between the lower (9 = 0) and upper (9 = 1) 
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stationary states (it is assumed that the surface tension decreases linearly with growing 
concentration). 

If the front motion is sufficiently slow, so that the flow field remains quasistationary, 
Eq. (|]) applies, which gives a closed equation of the front displacement. This equation 
is nonlocal so that the motion is dependent on the instantaneous distribution of domains 
occupied by the alternative states in the entire region. 

The thickness of the front region where the change of the surfactant concentration occurs 



is determined by the characteristic diffusional length y ID jr. On distances far exceeding this 
scale, the concentration distribution can be considered as stepwise. Using the discontinuous 
function #(x) that assumes the two alternative values 6 = and 9 = 1 in the domains 
separated by a boundary T (presumed smooth but not necessarily simply connected) brings 
Eq. @ to the form 

, x M f n , i > 

c x =-—f ^ds, 10 
Air Jr r(s) 

where the contour is parametrized by the arc length s, and r(s) = |x — £(s)| is the distance 
from the reference point to a point on the contour; by convention, the normal is directed 
towards the domain occupied by the upper state. 

Eq. ([ToD is quite adequate for computing the flow velocity far from the fronts but cannot 
be applied to the problem of front dynamics, since the integral diverges on the front location. 
To compute the velocity induced within the front region, one has to remove a short arc 
segment from the integral Eq. (|10"D, and apply on the rest of the contour Eq. @ containing 
a smooth concentration profile. On the removed segment, the finite diffusional length has 



taken into account, and both parts of the contour matched as in Ref fl4| . 

We presume that the local curvature radius of the front is of the same order of magnitude 
as the characteristic macroscopic scale L, and far exceeds the diffusional length. Then the 
concentration gradient is directed normally to the front, and the contour segments are almost 
rectilinear on distances s ^ L. In accordance to a common procedure of the theory of front 
dynamics H,|K|, we transform to a local coordinate frame comoving with the front, and take 



a certain contour, say, 6 = q as the origin of the normal coordinate £. The origin of the 
tangential coordinate rj is taken at the foot of the normal drawn from the reference point. 
Consider a point with a normal coordinate £ -C L, and cut from the contour integral in 
Eq. (|I0D a small arc of the length 25, where ^ < (5 < L. On the remaining part of the 
contour, £ can be neglected; thus, the velocity induced at the reference point by the "far" 
contour is determined by Eq. fllQ|) . The normal velocity computed in this way diverges 
logarithmically at 5 — > 0, and can be expressed as 

where C is a constant dependent on the shape of the contour in the far region. In the near 
region, we introduce a stretched coordinate T] = s/£, and compute the inner integral as 

M i- f°° nl/t- mj/ f S/C d V 



c«(£) = -- lim / e'ii-QdQ 

Z7T 5/^^oo J-oo JO 



M 

2^ 



/oo 
InlCl^e-CK • (12) 
-oo 



The auxilliary value 5 falls out when both integrals are added up. The resulting finite 
expression for the flow velocity in the front region is 

c(0 = -^ln^y; ln/3(0 = /jn|e-CI^'(CK. (13) 

Since #'(C) — » at ( — > oo, /3(£) x £ at £ ^> 1. The function /3(£) computed numerically 
using the analytical solution 0'(£) = |sech 2 | corresponding to Eq. (§) is plotted in Fig. |l|. 

In view of the nonlocal character of the induced motion, only a symmetric configuration 
may be stationary. Consider a spot of the radius L occupied by the lower state on the 
infinite interface occupied by the upper state. We presume that the spot radius is large 



when measured on the diffusional scale, so that the dimensionless radius I = L/ylDr 1. 
The induced surface flow velocity is computed using Eq. (|10D as 

. . Mr cos ada 

c{p) = -7T- 



2tt Jo 



sfl + (p/l) 2 - 2(p/l) cos a 
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at p > I, 



(14) 



at p < I. 

where p, a are polar coordinates, and K, E are complete elliptic integrals. The radial flow 
velocity vanishes at r = and decreases ocl/ratr^oo (Fig. 
As expected, there is a divergence in the vicinity of the front: 



M 
'2tt 



In 



e 2 e 



at £ 0, 



(15) 



where £ = p — I, and Z ^> £ = 0(1). The divergence is resolved using the asymptotic 
procedure described above. The resulting finite expression for the flow velocity at distances 
£ = 0(1) from the spot circumference coincides with Eq. ([13]) where C = 81 /e 2 . In the 
outer limit x — > oo, Eq. ( [13]) coincides with the inner limit ([15]) of the outer solution. The 
uniformly valid composite expansion can be obtained by adding up Eqs. (|T3"D , (p]) and 
extracting their common limit (|i~5]) : 



c(p) 



At 



Af 



(16) 



For the purpose of semi- quantitative estimation, one can assume that the front propaga- 
tion is mostly affected by the maximum advection velocity in the front region, which is ob- 
served at the location where the concentration gradient is also at its maximum. Using the nu- 
merical value (3(0) = 0.882, the maximum velocity is computed as c m = —(M/2ir) ln(0.955/). 
Since the dependence on I is logarithmical, this velocity may be of 0(1) on the diffusional 
scale, and balance the intrinsic front velocity Co = q — \ when I S> 1 but In I is not exceedingly 
large. The radius of a stationary spot is computed then as 



L = 1.047v / 2^exp 71(29 l) 

M 



(17) 



Since the variable flow velocity is actually smaller than its maximal value, a somewhat 
larger radius might be actually necessary to induce the flow counterbalancing the speed 



of the chemical front; thus, the above expression gives a lower estimate of the spot size; 
the correction would only amount, however, to modifying the numerical coefficient. The 
spot should be stable to homogeneous perturbations, since an increase of the radius would 
cause an increased inward flow restoring the stationary state. The curvature correction to 
the intrinsic speed of the chemical front is of 0(l~ r ), and can be neglected. I checked the 
stability using the method of boundary perturbations similar to that of Ref. HlQfl . A reverse 



set-up, with a spot occupied by the upper state (with a lower surface tension) is treated 
formally in the same way, but, as the flow in this case is outward, the stationary state is 
evidently unstable. A similar instability was observed in experiment with radial spreading 
of a surfactant [plj . 



A rational expansion can be carried out for the case when the Marangoni number is small, 
and the two alternative states are close to the Maxwell construction. The most interesting 
situation arises when front curvature effects are also of the same order of magnitude. Thus, 
we can take Z _1 as the small parameter of the expansion, and set q = ^ + qi/l, M = M\jl. 
Expanding 9 = 6q + l^ l 9\ + . . . , c = l~ l c\ + . . ., yields in the zero order 6$ defined by Eq. ([|). 
The stationary first order equation valid in the vicinity of the boundary of a circular spot 
with the radius I is 

\0{{p) + /'(Wi + ¥'o(p) - <?A(1 - ) ~ Mo)' = 0. (18) 

The homogeneous part of this equation has a zero eigenvalue corresponding to the trans- 
lational symmetry of the zero-order equation; the corresponding eigenfunction is 0' (p) = 
|sech 2 |; £ = p — I. The solvability condition of Eq. ( fL8|) requires that its inhomogeneous 
part be orthogonal to this eigenfunction. This yields a relation between q\ and M\ that has 
to be satisfied at a stationary front. Since the eigenfunction falls off sharply at £ ^> 1, the 
asymptotic expression for the velocity valid in the front region, Eq. ([13]) with M replaced 
by M\ can be used for c\. Computing the integrals, we express the solvability condition as 

1 Mi I 

2-" + 27 ln ir - (19) 



where 

(t)0o(t)d£ / In \£-C\e (QdC= 1-344. (20) 

-OO J—oo 

The first term in Eq. ([[£]) describes the curvature effect. If it is neglected, Eq. ( pf ) is 
recovered, but with a larger numerical coefficient exp(Z ) = 1,969. 

If there is a number of spots removed at distances far exceeding the diffusional length, the 
induced flow is additive, since the Stokes problem is linear. Separate spots are attracted one 
to the other by the induced inward flow; clearly, the system shows a tendency to aggregation 
but coalescing spots should srink again to the same stable radius at which the inward flow 
is compensated by the chemical propagation speed. Thus, Marangoni flow serves as an 
efficient mechanism preventing the spread of the state with a higher surface tension when it 
is dynamically favorable. 
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FIG. 2. The interfacial radial velocity in and around a solitary spot. 
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